* Do file to replicate the Tables in:
* Leon and Miguel (2016) Risky Transportation choices and the Value of Statistical Life," AEJ: Applied.

clear all
set more off
cap log close

* Set the directory where the data files are HERE
*do 	"/Users/Gianmarco/Dropbox/research projects/Transportation choices and the value of life/2012/Do/Direct_loc.do"
do 	"/Users/Gianmarco 1/Dropbox/research projects/Transportation choices and the value of life/2012/Do/Direct_loc-1.do"


* Table 1: Descriptive about the modes of transport
	* Use data in wide format (one observation per passenger-trip)
	use		"$Transp_data/Transp-wide-final-Replication.dta", clear

	* Risk
		* Risk of being in an accident and dying
		sum risk_obs1* [aw=weight] if sample==1
		* Risk of being in an accident
		sum risk_obs2* [aw=weight] if sample==1
	* Perceived qualities of the different modes
   		* Sufix 1= Helicopter ; Sufix 2=Water Taxi ; Sufix 3=Ferry ; Sufix 4 = Hovercraft
		sum seats* noise* crowdness* convenient_loc* clientele* if sample==1		

* Table 3: Respondent descriptive statistics

* Panel A: Transport Choices
	* Use dataset in long format (passengerXtrip)
	use		"$Transp_data/Transp-Regressions.dta", clear
	tab mode, g(mode_)
	* Mode 1=Helicopter ; Mode 2=Water Taxi ; Mode 3=Ferry ; Mode 4= Hovercraft
	
		* keep one obs per passenger trip
	keep if alt==1
	
	sum 		mode_1 mode_2 mode_3 mode_4 [w=weight] if sample_af==1
	sum 		mode_1 mode_2 mode_3 mode_4 [w=weight] if sample_noaf==1
	sum 		mode_1 mode_2 mode_3 mode_4 [w=weight] if sample_all==1
	
* Panel B: Respondent Characteristics and Attitudes
	use		"$Transp_data/Transp-wide-final-Replication.dta", clear

	sum 		gender age educ_low educ_high affected children swim nationality_1 nationality_2 nationality_3  wage_hr_ppp wage_imput_ppp life_exp fatalism [w=weight] if nationality<3 & sample==1/*African*/
	sum 		gender age educ_low educ_high affected children swim nationality_1 nationality_2 nationality_3  wage_hr_ppp wage_imput_ppp life_exp fatalism [w=weight] if  nationality==3 & sample==1/*Non-African*/
	sum 		gender age educ_low educ_high affected children swim nationality_1 nationality_2 nationality_3  wage_hr_ppp wage_imput_ppp life_exp fatalism [w=weight] if sample==1/*African*/

* Table 4: Transportation Choices and the Value of a Statistical Life � Mixed logit estimates
	*See MATLAB replication files:
	* mixed_logit_RestTriang_20160620_T4_1
	* mixed_logit_RestTriang_20160620_T4_2
	
* Table 5: Value of a Statistical Life estimates in different subsamples and with alternative specifications 

	* For Mixed Logit estimates (Rows 1-10), see MATLAB replication files:
	* mixed_logit_RestTriang_20160620_T5_1
	* mixed_logit_RestTriang_20160620_T5_2

	* Conditional logit estimates: Rows 11-13
	use		"$Transp_data/Transp-Regressions.dta", clear

	global attributes="seats noise crowdness convloc clientele"

	clogit	alt  risk_obs transp_cost $attributes [iw= weight] if sample_all==1, group(choice_situation) cl(serial) technique(bfgs)
			outreg	using "CLogit", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\N_clust\ll) summtitle("Observations"\"Num. of decision makers"\"Log-Likelihood") replace nodisp
			codebook choice_situation if e(sample)
			nlcom -(_b[risk_obs])/(_b[transp_cost])

	clogit	alt  risk_obs transp_cost $attributes [iw= weight] if sample_af==1, group(choice_situation) cl(serial) technique(bfgs)
			outreg	using "CLogit", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\N_clust\ll) summtitle("Observations"\"Num. of decision makers"\"Log-Likelihood") merge nodisp
			codebook choice_situation if e(sample)
			nlcom -(_b[risk_obs])/(_b[transp_cost])

	clogit	alt  risk_obs transp_cost $attributes [iw= weight] if sample_noaf==1, group(choice_situation) cl(serial) technique(bfgs)
			outreg	using "CLogit", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\N_clust\ll) summtitle("Observations"\"Num. of decision makers"\"Log-Likelihood") merge nodisp
			codebook choice_situation if e(sample)
			nlcom -(_b[risk_obs])/(_b[transp_cost])

* Table 6: Determinants of the VSL across individuals

	**** Dataset that incorporates the individual level estimates of the VSL from the Mixed Logit regressions
	use "$Transp_data/Trans_VSLEstimates.dta", clear

	gen		l_vsl=log(vsl)
	gen		l_wage=log(wage_imput_ppp)
	g		african=(nationality<=2)

	bootstrap, reps(1000): reg 	l_vsl  african  if l_wage~=., r
			outreg	using "VSLDeterm-2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\r2) summtitle("Observations"\"R-sq") replace nodisp
	bootstrap, reps(1000): reg 	l_vsl  african life_exp if l_wage~=., r
			outreg	using "VSLDeterm-2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\r2) summtitle("Observations"\"R-sq") merge nodisp
	bootstrap, reps(1000): reg 	l_vsl  african fatalism if l_wage~=., r
			outreg	using "VSLDeterm-2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\r2) summtitle("Observations"\"R-sq") merge nodisp
	bootstrap, reps(1000): reg 	l_vsl	l_wage, r
			outreg	using "VSLDeterm-2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\r2) summtitle("Observations"\"R-sq") merge nodisp
	bootstrap, reps(1000): reg 	l_vsl  african l_wage, r
			outreg	using "VSLDeterm-2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\r2) summtitle("Observations"\"R-sq") merge nodisp
	bootstrap, reps(1000): reg 	l_vsl  african life_exp fatalism l_wage , r
			outreg	using "VSLDeterm-2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\r2) summtitle("Observations"\"R-sq") merge

* Table A.1: Descriptive Statistics, by mode of transportation
	use		"$Transp_data/Transp-wide-final-Replication.dta", clear
	
	sum 		safer cheaper faster comfortable gender age educ_low educ_high affected children swim nationality_1 nationality_2 nationality_3 wage_hr_ppp  wage_imput_ppp life_exp  fatalism [w=weight]  if mode==3 & sample_all==1 /* Ferry*/
	sum 		safer cheaper faster comfortable gender age educ_low educ_high affected children swim nationality_1 nationality_2 nationality_3 wage_hr_ppp  wage_imput_ppp life_exp  fatalism [w=weight]  if mode==2 & sample_all==1 /* Water Taxi */
	sum 		safer cheaper faster comfortable gender age educ_low educ_high affected children swim nationality_1 nationality_2 nationality_3 wage_hr_ppp  wage_imput_ppp life_exp  fatalism [w=weight]  if mode==4 & sample_all==1 /* Hovercraft */

* Table A.2: Transportation Choices and the Value of a Statistical Life: Robustness Checks
	* See MATLAB Replication files:
	* mixed_logit_RestTriang_20160620_T5_1

* Table A.3: Transportation Choices and the Value of a Statistical Life (conditional logit)
	use		"$Transp_data/Transp-Regressions.dta", clear 

	global attributes="seats noise crowdness convloc clientele"
	global africa="risk_obs_africa transp_cost_africa"

	clogit	alt  risk_obs transp_cost [iw= weight] if sample_af==1, group(choice_situation) cl(serial) technique(bfgs)
			outreg	using "CLogit2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\N_clust\ll) summtitle("Observations"\"Num. of decision makers"\"Log-Likelihood") replace nodisp
			codebook choice_situation if e(sample)
			nlcom -(_b[risk_obs])/(_b[transp_cost])
	clogit	alt  risk_obs transp_cost $attributes [iw= weight] if sample_af==1, group(choice_situation) cl(serial) technique(bfgs)
			outreg	using "CLogit2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\N_clust\ll) summtitle("Observations"\"Num. of decision makers"\"Log-Likelihood") merge nodisp
			codebook choice_situation if e(sample)
			nlcom -(_b[risk_obs])/(_b[transp_cost])

	clogit	alt  risk_obs transp_cost [iw= weight] if sample_noaf==1, group(choice_situation) cl(serial) technique(bfgs)
			outreg	using "CLogit2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\N_clust\ll) summtitle("Observations"\"Num. of decision makers"\"Log-Likelihood") merge nodisp
			codebook choice_situation if e(sample)
			nlcom -(_b[risk_obs])/(_b[transp_cost])
	clogit	alt  risk_obs transp_cost $attributes [iw= weight] if sample_noaf==1, group(choice_situation) cl(serial) technique(bfgs)
			outreg	using "CLogit2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\N_clust\ll) summtitle("Observations"\"Num. of decision makers"\"Log-Likelihood") merge nodisp
			codebook choice_situation if e(sample)
			nlcom -(_b[risk_obs])/(_b[transp_cost])

	clogit	alt  risk_obs transp_cost [iw= weight] if sample_all==1, group(choice_situation) cl(serial) technique(bfgs)
			outreg	using "CLogit2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\N_clust\ll) summtitle("Observations"\"Num. of decision makers"\"Log-Likelihood") merge nodisp
			codebook choice_situation if e(sample)
			nlcom -(_b[risk_obs])/(_b[transp_cost])
	clogit	alt  risk_obs transp_cost $attributes [iw= weight] if sample_all==1, group(choice_situation) cl(serial) technique(bfgs)
			outreg	using "CLogit2", starlevels(10 5 1) bdec(3) varlabels  se summstat(N\N_clust\ll) summtitle("Observations"\"Num. of decision makers"\"Log-Likelihood") merge
			codebook choice_situation if e(sample)
			nlcom -(_b[risk_obs])/(_b[transp_cost])

* Table A.4: Transportation Choices and the Value of a Statistical Life � Mixed logit estimates
		* See MATLAB Replication files:
		* mixed_logit_RestTriang_20160620_TA4


**** FIGURES

* Figure 1: Map of the Study Setting, Lungi International Airport and Freetown, Sierra Leone
	* Google Maps and own elaboration
	
 
* Figure 2: Operation and pricing for the modes of transport 
	* Administrative data

* Figure 3: Transport choice as a function of wages (w) and value of life (VSL)
	* Illustration
	
* Figure 4: Perceived Transportation Risk Rankings 
	
	*(observation per passenger, always for the rainy season)
	use		"$Transp_data/Transp-wide-final-Replication.dta", clear
	
	* Sufix 1= Helicopter ; Sufix 2=Water Taxi ; Sufix 3=Ferry ; Sufix 4 = Hovercraft
	tab safest1_rain [aw=weight] if sample==1
	tab safest2_rain [aw=weight] if sample==1
	tab safest3_rain [aw=weight] if sample==1
	tab safest4_rain [aw=weight] if sample==1

* Figure 5: Distribution of individual VSL estimates, mixed logit estimates with restricted triangular distributions

	**** Dataset that incorporates the individual level estimates of the VSL from the Mixed Logit regressions
	use "$Transp_data/Trans_VSLEstimates.dta", clear

	* PANEL A: Density for the Full Sample: 
	sum vsl, d
	kdensity vsl if vsl<=r(p99), g(x y) nogr
	label var y "Full Sample"
	twoway (line y x  if x<=1200, lcolor(black) lwidth(thick) lpattern(solid)), legend(size(small) ring(0) cols(1) position(3)) ytitle(Density) xtitle(Value of Statistical Life (in '000 US$ PPP)) title(Panel A: Full Sample, size(medium) margin(vsmall)) graphregion(fcolor(white))  saving(FullSample, replace)
	 
	* Panel B: African Non African                                             
	kdensity vsl if nationality<=2, g(y2) at(x) nogr
	kdensity vsl if nationality==3, g(y3) at(x) nogr
				label var y2 "Africans"	
				label var y3 "Non-Africans"	

	twoway (line y2 x if x<=1200, lcolor(black) lpattern(solid)) (line y3 x if x<=1200, lcolor(black) lpattern(dash)), legend(size(small) ring(0) cols(1) position(3)) ytitle(Density) xtitle(Value of Statistical Life (in '000 US$ PPP)) title(Panel B: African and Non-Africans, size(medium) margin(vsmall)) graphregion(fcolor(white))  saving(Africans, replace)
	drop	y2 y3

	* Panel C: Full sample, by wage quartiles

	sum wage_imput_ppp, d
	kdensity vsl if wage_imput_ppp<r(p25), g(y1) at(x) nogr 
	sum wage_imput_ppp, d
	kdensity vsl if wage_imput_ppp>=r(p75), g(y3) at(x) nogr
				label var y1 "Poorest Quartile (<US$13.5)"	
				label var y3 "Richest Quartile (>US$48)"	
	twoway (line y1 x if x<=1200, lcolor(black) lpattern(solid) legend(size(vsmall)) sort) (line y3 x  if x<=1200, lcolor(black) lpattern(dash) legend(size(vsmall) ring(0) cols(1) position(3)) sort), legend(size(small)) ytitle(Density) xtitle(Value of Statistical Life (in '000 US$ PPP)) title("Panel C: Full Sample, by Wage Quartiles", size(medium) margin(vsmall)) graphregion(fcolor(white))  saving(AllIncome, replace)
	drop	y1 y3

	* Panel D: Full Sample, By fatalism
	sum fatalism, d
	kdensity vsl if fatalism<=4, g(y1) at(x) nogr 
	sum fatalism, d
	kdensity vsl if fatalism>=7, g(y2) at(x) nogr
				label var y1 "Non-Fatalistic (<=5)"	
				label var y2 "Fatalistic (>=6)"	
	twoway (line y1 x if x<=1200, lcolor(black) lpattern(solid) sort) (line y2 x if x<=1200, lcolor(black) lpattern(dash) sort), legend(size(small) ring(0) cols(1) position(3)) ytitle(Density) xtitle(Value of Statistical Life (in '000 US$ PPP)) title("Panel D: Full Sample, by Fatalism", size(medium) margin(vsmall)) graphregion(fcolor(white))  saving(AllFatalism, replace)
	drop	y1 y2

	graph combine 	FullSample.gph Africans.gph AllIncome.gph AllFatalism.gph, ycommon xcommon ysize(8) xsize(12) row(2) col(2) graphregion(fcolor(white)) imargin(small) saving(Figure5_VSL, replace)

* Figure A.1: Timing of the trips available in the dataset
	use		"$Transp_data/Transp-Regressions.dta", clear 
	* keep one obs per passenger trip
	keep if alt==1
	g trim=1 if month>=1 & month<=3
		replace trim=2 if  month>=4 & month<=6
		replace trim=3 if  month>=7 & month<=9
		replace trim=4 if  month>=10 & month<=12
	
	g 	year_trim=year*10+trim
	
	tab	year_trim
	
* Figure A.2: Observable Differences between African Travelers and non-African Travelers
	use		"$Transp_data/Transp-wide-final-Replication.dta", clear

	* Panel A: Overall Quality: Good or Excellent                          
	foreach num of numlist 1/4 {
			gen good`num'=(overall_qual`num'>3)
		bys africa_born: tab good`num' [aw=weight]
		}
		
	* Panel B: Hourly Wage
	kdensity wage_imput_ppp if wage_imput_ppp<250, g(x y) nogr
	kdensity wage_imput_ppp if wage_imput_ppp<250 & nationality==1, g(y1) at(x) nogr 
	kdensity wage_imput_ppp if wage_imput_ppp<250 & africa_born==1, g(y2) at(x) nogr
	kdensity wage_imput_ppp if wage_imput_ppp<250 & africa_born==0, g(y3) at(x) nogr
				label var y1 "Sierra Leone"	
				label var y2 "Africans"	
				label var y3 "Non-Africans"	
	twoway (line y2 x, lcolor(black) lpattern(solid)) (line y3 x, lcolor(black) lpattern(dash)), legend(size(small) ring(0) cols(1) position(3)) ytitle(Density) xtitle(Hourly Wage (Imputed) - USD PPP) saving(Wages, replace)
	drop	x y y1 y2 y3

	* Panel C: Self Reported Life Expectancy
	kdensity life_exp, g(x y) nogr
	kdensity life_exp if nationality==1, g(y1) at(x) nogr 
	kdensity life_exp if africa_born==1, g(y2) at(x) nogr
	kdensity life_exp if africa_born==0, g(y3) at(x) nogr
				label var y1 "Sierra Leone"	
				label var y2 "Africans"	
				label var y3 "Non-Africans"	
	twoway (line y2 x, lcolor(black) lpattern(solid)) (line y3 x, lcolor(black) lpattern(dash)), xlabel(5(5)70) xsc(r(5 70))  legend(size(small) ring(0) cols(1) position(9)) ytitle(Density) xtitle(Self Reported Remaining Life Expectancy) saving(LifeExpect, replace)

	drop	x y y1 y2 y3

	* Panel D: Fatalism
	bys africa_born: tab fatalism [aw=weight]

* Figure A.3: Example and definition of a Restricted Triangular Distribution
	* Illustration




		


		
